Uniform convergence guarantees for the deep Ritz method for nonlinear problems

We provide convergence guarantees for the Deep Ritz Method for abstract variational energies. Our results cover nonlinear variational problems such as the p-Laplace equation or the Modica–Mortola energy with essential or natural boundary conditions. Under additional assumptions, we show that the convergence is uniform across bounded families of right-hand sides.


Introduction
The idea of the Deep Ritz Method is to use variational energies as an objective function for neural network training to obtain a finite-dimensional optimization problem that allows solving the underlying partial differential equation approximately. The idea of deriving a finite-dimensional optimization problem from variational energies dates back to Ritz [28], was widely popularized in the context of finite element methods (see, e.g., Braess [4]), and was recently revived by E and Yu [13] using deep neural networks. In the following, we give a more thorough introduction to the Deep Ritz Method. Let ⊆ R d be a bounded domain and consider the variational energy corresponding to the Lagrangian L and a force f , namely defined on a suitable function space X, usually a Sobolev space W 1,p ( ). One is typically interested in minimizers of E on subsets U ⊆ X where U encodes further physical constraints, such as boundary conditions. Here, we consider either unconstrained problems or zero Dirichlet boundary conditions and use the notation U = X 0 for the latter case. In other words, for zero boundary conditions, one aims to find To solve such a minimization problem numerically, the idea dating back to Ritz [28] is to use a parametric ansatz class and to consider the finite-dimensional minimization problem of finding which can be approached by different strategies, depending on the class A. For instance, if A is chosen to be a finite element ansatz space or polynomials and the structure of E is simple enough, one uses optimality conditions to solve this problem.
In this manuscript, we focus on ansatz classes that are given through (deep) neural networks. When choosing such ansatz functions, the method is known as the Deep Ritz Method and was recently proposed by E and Yu [13]. Neural network type ansatz functions possess a parametric form as in (3), however, it is difficult to impose zero boundary conditions on the ansatz class A. To circumvent this problem, one can use a penalty approach, relaxing the energy to the full space, but penalizing the violation of zero boundary conditions, to include these. This means that for a penalization parameter λ > 0 one aims to find The idea of using neural networks for the approximate solution of PDEs can be traced back at least to the works of Lee and Kang [21], Dissanayake and Phan-Thien [10], Takeuchi and Kosugi [32], Lagaris et al. [20]. Since the recent successful application of neural network based methods to stationary and instationary PDEs by E et al. [12], E and Yu [13], Sirignano and Spiliopoulos [30], there is an ever growing body of theoretical works contributing to the understanding of these approaches. For a collection of the different methods, we refer to the overview articles by Beck et al. [3], Han et al. [15].
The error in the Deep Ritz Method, which decomposes into an approximation, optimization, and generalization terms, has been studied by Luo and Yang [25], Xu [34], Duan et al. [11], Hong et al. [17], Jiao et al. [18], Lu et al. [23], Lu et al. [24], Müller and Zeinhofer [26]. However, those works either consider non-essential boundary conditions or they require a term with a positive potential, apart from Müller and Zeinhofer [26]. This excludes the prototypical Poisson equation, which was originally treated by the Deep Ritz Method by E and Yu [13]. More importantly, those works only study linear problems, which excludes many important applications.
In this work, we thus study the convergence of the Deep Ritz Method when a sequence of growing ansatz classes (A n ) n∈N , given through parameter sets n and a penalization of growing strength (λ n ) n∈N with λ n ∞, is used in the optimization problem (4) with more modest assumptions on L, f , and .
Denote a sequence of (almost) minimizing parameters of problem (4) with parameter set n and penalization λ n by θ n . We then see that under mild assumptions on (A n ) n∈N and E, the sequence (u θ n ) n∈N of (almost) minimizers converges weakly in X to the solution of the continuous problem, see Theorem 7 in Sect. 3. We then strengthen this result in Sect. 4 where we show that the aforementioned convergence is uniform across certain bounded families of right-hand sides f , see Theorem 12. This means that a fixed number of degrees of freedom in the ansatz class can be used independently of the right-hand side to achieve a given accuracy. Alternatively, given a discretization of the space of right-hand sides, one may discretize the solution operator that maps f to the minimizer u of (2) and still obtain a convergence guarantee (although this is not necessarily a viable numerical approach).
To the best of our knowledge, our results currently comprise the only convergence guarantees for the Deep Ritz Method for nonlinear problems. However, since we prove these results using -convergence methods, no rates of convergence are obtained -as mentioned above, for linear elliptic equations some error decay estimates are known. Our results also do not provide insight into the finite-dimensional optimization problem (4) which is a challenging problem in its own right, see, for instance, Wang et al. [33], Courte and Zeinhofer [8]. However, they guarantee that given one is able to solve (4) to a reasonable accuracy, one is approaching the solution of the continuous problem (2).
Our results are formulated for neural network type ansatz functions due to the current interest in using these in numerical simulations, yet other choices are possible. For instance, our results do apply directly to finite element functions.
The remainder of this work is organized as follows. Section 2 discusses some preliminaries and the used notation. The main results, namely -convergence and uniformity of convergence are provided in Sects. 3 and 4, respectively. Finally, in Sect. 5 we discuss how the p-Laplace and a phase field model fit into our general framework.

Notation and preliminaries
We fix our notation and present the tools that our analysis relies on.

Notation of Sobolev spaces and Friedrich's inequality
We denote the space of functions on ⊆ R d that are integrable in the pth power by L p ( ), where we assume that p ∈ [1, ∞). Endowed with u p L p ( ) := |u| p dx, this is a Banach space, i.e., a complete normed space. If u is a multivariate function with values in R m , we interpret |·| as the Euclidean norm. We denote the subspace of L p ( ) of functions with weak derivatives up to order k in L p ( ) by W k,p ( ), which is a Banach space with the norm This space is called a Sobolev space and we denote its dual space, i.e., the space consisting of all bounded and linear functionals on W k,p ( ) by W k,p ( ) * . The closure of all compactly supported smooth functions C ∞ c ( ) in W k,p ( ) is denoted by W k,p 0 ( ). It is well known that if has a Lipschitz continuous boundary the operator that restricts a Lipschitz continuous function on to the boundary admits a linear and bounded extension tr : W 1,p ( ) → L p (∂ ). This operator is called the trace operator and its kernel is precisely W 1,p 0 ( ). Further, we write u L p (∂ ) whenever we mean tr(u) L p (∂ ) . In the following we mostly work with the case p = 2 and write H k (0) ( ) instead of W k,2 (0) ( ). In order to study the boundary penalty method, we use the Friedrich inequality which states that the L p ( ) norm of a function can be estimated by the norm of its gradient and boundary values. We refer to Gräser [14] for a proof.

Neural networks
Here we introduce our notation for the functions represented by a feedforward neural network. Consider natural numbers d, m, L, N 0 , . . . , N L ∈ N and let The neural network function with parameters θ and with respect to some activation function ρ : R → R is the function The set of all neural network functions of a certain architecture is given by {u ρ θ | θ ∈ }, where collects all parameters of the above form with respect to fixed natural numbers d, m, L, N 0 , . . . , N L . If we have f = u ρ θ for some θ ∈ we say the function f can be realized by the neural network F ρ . Note that we often drop the superscript ρ if it is clear from the context.
A particular activation function often used in practice and relevant for our results is the rectified linear unit or ReLU activation function, which is defined via x → max{0, x}. Arora et al. [2] showed that the class of ReLU networks coincides with the class of continuous and piecewise linear functions. In particular, they are weakly differentiable. Since piecewise linear functions are dense in H 1 0 ( ), we obtain the following universal approximation result which we prove in detail in the appendix.

Theorem 2 (Universal approximation with zero boundary values) Consider an open set
⊆ R d and fix a function u ∈ W 1,p 0 ( ) with p ∈ [1, ∞). Then for all ε > 0 there exists u ε ∈ W 1,p 0 ( ) that can be realized by an ReLU network of depth log 2 (d + 1) + 1 such that To the best of our knowledge, this is the only available universal approximation result where the approximating neural network functions are guaranteed to have zero boundary values. This relies on the special properties of the ReLU activation function and it is unclear for which classes of activation functions universal approximation with zero boundary values hold.

Gamma convergence
We recall the definition of -convergence with respect to the weak topology of reflexive Banach spaces. For further reading, we point the reader towards Dal Maso [9]. Definition 3 ( -convergence) Let X be a reflexive Banach space as well as F n , F : X → (-∞, ∞]. Then (F n ) n∈N is said to be -convergent to F if the following two properties are satisfied: is bounded in X (or equivalently, relatively compact with respect to the weak topology) for all r ∈ R. We say that a sequence (x n ) n∈N are quasiminimizers of the functionals (F n ) n∈N if we have where δ n → 0.
We need the following property of -convergent sequences. We want to emphasize the fact that there are no requirements regarding the continuity of any of the functionals and that the functionals (F n ) n∈N are not assumed to admit minimizers.

Theorem 4 (Convergence of quasiminimizers) Let X be a reflexive Banach space and
(F n ) n∈N be an equicoercive sequence of functionals that -converges to F. Then, any sequence (x n ) n∈N of quasiminimizers of (F n ) n∈N is relatively compact with respect to the weak topology of X and every weak accumulation point of (x n ) n∈N is a global minimizer of F. Consequently, if F possesses a unique minimizer x, then (x n ) n∈N converges weakly to x.

Abstract -convergence result for the deep Ritz method
For the abstract results, we work with an abstract energy E : X → R, instead of an integral functional of the form (1). This reduces technicalities in the proofs and separates abstract functional-analytic considerations from applications.

Setting 5
Let (X, · X ) and (B, · B ) be reflexive Banach spaces and γ ∈ L(X, B) be a continuous linear map. We set X 0 to be the kernel of γ , i.e., X 0 = γ -1 ({0}). Let ρ : R → R be some activation function and denote by ( n ) n∈N a sequence of neural network parameters.

We assume that any function represented by such a neural network is a member of X and define
Here, x θ denotes the function represented by the neural network with the parameters θ . Let E : X → (-∞, ∞] be a functional and (λ n ) n∈N a sequence of real numbers with λ n → ∞. Furthermore, let p ∈ (1, ∞) and f ∈ X * be fixed and define the functional Then assume the following holds: The functional E is bounded from below, weakly lower semicontinuous with respect to the weak topology of (X, · X ) and continuous with respect to the norm topology of (X, · X ). (A3) The sequence (F f n ) n∈N is equicoercive with respect to the norm · X .
Remark 6 We discuss Assumptions (A1) to (A3) in view of their applicability to concrete problems.
(i) In applications, (X, · X ) will usually be a Sobolev space with its natural norm, the space B contains boundary values of functions in X and the operator γ is a boundary value operator, e.g., the trace map. However, if the energy E is coercive on all of X, i.e., without adding boundary terms to it, we might choose γ = 0 and obtain X 0 = X. This is the case for non-essential boundary value problems. (ii) Assumption (A1) compensates that, in general, we cannot penalize with arbitrary strength. However, if we can approximate any member of X 0 by a sequence x θ n ∈ A n ∩ X 0 then any divergent sequence (λ n ) n∈N can be chosen. This is, for example, the case for the ReLU activation function and the space X 0 = H 1 0 ( ). More precisely, we can choose A n to be the class of functions expressed by a (fully connected) ReLU network of depth log 2 (d + 1) + 1 and width n, see Theorem 2.
Theorem 7 ( -convergence) Assume we are in Setting 5. Then the sequence (F f n ) n∈N of functionals -converges towards F f . In particular, if (δ n ) n∈N is a sequence of nonnegative real numbers converging to zero, any sequence of δ n -quasiminimizers of F f n is bounded and all its weak accumulation points are minimizers of F f . If, additionally, F f possesses a unique minimizer x f ∈ X 0 , any sequence of δ n -quasiminimizers converges to x f in the weak topology of X.
Proof We begin with the limes inferior inequality. Let x n x in X and assume that x / ∈ X 0 . Then f (x n ) converges to f (x) as real numbers and γ (x n ) converges weakly to γ (x) = 0 in B. Combining this with the weak lower semicontinuity of · p B , we get, using the boundedness from below, that Now let x ∈ X 0 . Then by the weak lower semicontinuity of E, we find Now let us have a look at the construction of the recovery sequence. For x / ∈ X 0 , we can choose the constant sequence and estimate A sufficient criterion for equicoercivity of the sequence (F f n ) n∈N from Assumption (A3) in terms of the functional E is given by the following lemma.
Lemma 8 (Criterion for equicoercivity) Assume we are in Setting 5. If there is a constant c > 0 such that it holds for all x ∈ X that Proof It suffices to show that the sequence So let r ∈ R be given and assume that r ≥ G f n (x). We estimate, assuming without loss of generality that λ n ≥ 1, As p > 1, a scaled version of Young's inequality clearly implies a bound on the set and hence the sequence (F f n ) n∈N is seen to be equicoercive.

Abstract uniform convergence result for the deep Ritz method
In this section we present an extension of Setting 5 that allows proving uniform convergence results over certain bounded families of right-hand sides.
Setting 9 Assume we are in Setting 5. Furthermore, let there be an additional norm |·| on X such that the dual space (X, |·|) * is reflexive. However, we do not require (X, |·|) to be complete. Then, let the following assumptions hold: (A4) The identity Id : (X, · X ) → (X, |·|) is completely continuous, i.e., maps weakly convergent sequences to strongly convergent ones. (A5) For every f ∈ X * , there is a unique minimizer x f ∈ X 0 of F f and the solution map is demicontinuous, i.e., maps strongly convergent sequences to weakly convergent ones.
Remark 10 As mentioned earlier, (X, · X ) is usually a Sobolev space with its natural norm. The norm |·| may then chosen to be an L p ( ) or W s,p ( ) norm, where s is strictly smaller than the differentiability order of X. In this case, Rellich's compactness theorem provides Assumption (A4).
Proof We begin by clarifying what we mean by S being defined on (X, |·|) * . Denote by i the inclusion map i : X 0 → X and consider By abusing notation, always when we refer to S as defined on (X, |·|) * we mean the above composition, i.e., Id •S • i * • Id * . Having explained this, it is clear that it suffices to show that Id * maps weakly convergent sequences to strongly convergent ones since i * is continuous, S demicontinuous, and Id strongly continuous. This, however, is a consequence of Schauder's theorem, see, for instance, Alt [1], which states that a linear map L ∈ L(X, Y ) between Banach spaces is compact if and only if L * ∈ L(Y * , X * ) is. Here, compact means that L maps bounded sets to relatively compact ones. Let X c denote the completion of (X, |·|). Then, using the reflexivity of (X, · X ) it is easily seen that Id : (X, · X ) → X c is compact. Finally, using that (X, |·|) * = X * c the desired compactness of Id * is established.
The following theorem is the main result of this section. It shows that the convergence of the Deep Ritz method is uniform on bounded sets in the space (X, |·|) * . The proof of the uniformity follows an idea from Cherednichenko et al. [7], where in a different setting a compactness result was used to amplify pointwise convergence to uniform convergence across bounded sets, compare to Theorem 4.1 and Corollary 4.2 in Cherednichenko et al. [7].
Theorem 12 (Uniform convergence of the Deep Ritz Method) Assume that we are in Setting 9 and let δ n 0 be a sequence of real numbers. For f ∈ X * , we set which is the approximate solution set corresponding to f and δ n . Furthermore, denote the unique minimizer of F f in X 0 by x f and fix R > 0. Then we have In the definition of this supremum, f is measured in the norm of the space (X, |·|) * . This means that f : (X, |·|) → R is continuous, which is a more restrictive requirement than the continuity with respect to · X . Also the computation of this norm takes place in the unit ball of (X, |·|), i.e., Before we prove Theorem 12, we need a -convergence result similar to Theorem 7. The only difference is that now also the right-hand side may vary along the sequence.

Proposition 13
Assume that we are in Setting 9, however, we do not need Assumption (A5) for this result. Let f n , f ∈ (X, |·|) * be such that f n f in the weak topology of the reflexive space (X, |·|) * . Then the sequence (F f n n ) n∈N of functionals -converges to F f in the weak topology of (X, · X ). Furthermore, the sequence (F f n n ) n∈N is equicoercive.
Proof The proof is almost identical to that of Theorem 7 but, since it is brief, we include it for the reader's convenience. We begin with the limes inferior inequality. Let x n x in X and x / ∈ X 0 . Then x n → x with respect to |·| which implies that f n (x n ) converges to f (x).
Using that γ (x n ) γ (x) in B, combined with the weak lower semicontinuity of · p B , we get Now let x ∈ X 0 . Then by the weak lower semicontinuity of E, we find Now let us have a look at the construction of the recovery sequence. For x / ∈ X 0 , we can choose the constant sequence and estimate As f n (X,|·|) * is bounded we find F f n n (x) → ∞ = F f (x). If x ∈ X 0 , we approximate it with a sequence (x n ) ⊆ X, according to Assumption (A1), such that x n ∈ A n and x n → x in · X and λ n γ (x n ) p B → 0. It follows that The equicoercivity was already assumed in (A3) so it does not need to be shown.

Proof of Theorem 12
We can choose (f n ) ⊆ (X, |·|) * and f n (X,|·|) * ≤ R and x f n n ∈ S n (f n ) such that Now it suffices to show that |x f n nx f n | converges to zero. Since (f n ) n∈N is bounded in (X, |·|) * and this space is reflexive, we can, without loss of generality, assume that f n f in (X, |·|) * . This implies by Lemma 11 that x f n → x f in (X, |·|). The -convergence result of the previous proposition yields x f n n x f in X and hence x f n n → x f with respect to |·| which concludes the proof.

Examples
We discuss different concrete examples that allow the application of our abstract results and focus on nonlinear problems. In particular, we consider a phase field model illustrating the basic -convergence result of Sect. 3 and the p-Laplacian as an example for the uniform results of Sect. 4.

General practical considerations
In practice, when solving the optimization problem (4), in order to obtain an approximate solution of the variational problem (2) there are a lot of choices to make. We give an overview over some of them here and report our specific choices in the individual examples.
Optimization One can use almost any kind of optimization algorithm for the approximation solution of the optimization problem (4), where gradient type algorithms and quasi-Newton methods are the most common choice. We use a combination of the Adam optimizer and L-BFGS. The former is a version of stochastic gradient descent with with adaptive moment estimation, see Kingma and Ba [19] and the latter is a quasi-Newton method, see Liu and Nocedal [22]. Typically, the optimization process begins with applying Adam until convergence slows down and the fast local convergence properties of Newton methods can be exploited through the application of the L-BFGS optimizer.
Quadrature In practice, one does not have access to the true gradient of the objective function. Hence, one usually uses estimates for the gradient as update directions. For example, E and Yu [13] used an online SGD estimator, i.e., a Monte Carlo approximation of the integral with fixed sample size. However, one can use any quadrature rule for the evaluation of the integrals in order to obtain approximations of the true gradient. We used a uniform grid for the discretization of the integral, i.e., the integral is approximated by the sum of the functions values at the grid points divided by the number of grid points. On the boundary of the domains, equispaced integration points are used and the quadrature rule is analogue to the one previously described.
We choose the number of integration points such that no further improved accuracy of the method is observed when increasing their amount. As this was computationally tractable without problems, no more elaborate integration routines were deemed necessary.
Activation functions The only requirements on the activation function present in Setting 5 are that the associated neural network functions belong to the considered Banach space A n ⊆ X as well as that Condition (A1) holds. The first one is usually of no concern as in practice X is often a space of given smoothness and hence for sufficently smooth activations A n ⊆ X is satisfied. Note that this is in particular the case X = H 1 ( ) and for the ReLU activation function. Further, in order to Condition (A1) to hold, it is necessary that the neural network type ansatz classes (A n ) have the universal approximation property in X 0 , i.e., that X 0 ⊆ n∈N A n . Note that this is the case for shallow 1 networks of increasing width and X = H k ( ) as long as the activation function is k times continuously differentiable and nonpolynomial (Pinkus [27]).
Penalization strength Condition (A1) couples the penalization strength to the norm of the (generalized) boundary values required for approximation of a general element x ∈ X. Consider, for example, the case that for any x ∈ X 0 there are x n ∈ X n such that x n → x and γ (x n ) B ≤ c(x)δ n for some δ n → 0. Then any choice of penalization strengths λ n → ∞ with λ n δ n → 0 satisfies Condition (A1). Let us first consider the case with inessential boundary values, which corresponds to γ ≡ 0 in the notation of Setting 5. Then, as argued above smooth and nonpolynomial activations together with arbitrarily strong penalization strengths are allowed.
Let us consider the case X = H 1 ( ) and γ = tr and B = L 2 ( ). For the ReLU activation function, Theorem 2 guarantees the existence of the u n → u with u n L 2 (∂ ) = 0, which allows for arbitrarily strong penalization. For other activation functions, which do not possess the universal approximation property with exact (generalized) zero boundary values, the proof of (A1) is more delicate and has to be established in specific cases.

A phase field model
Let ε > 0 be fixed, ⊆ R d a bounded Lipschitz domain and consider the following energy: where W : R → R is the nonlinear function given by The functional E constitutes a way to approximately describe phase separation and the parameter ε encodes the length-scale of the phase transition, see Cahn and Hilliard [6].
We describe now how the Setting 5 is applicable to fully connected neural network ansatz functions with tanh activation function. For the Banach spaces in Setting 5, we choose We choose γ ≡ 0, hence the choice of the space B and its norm is irrelevant. The choice of γ ≡ 0 corresponds to the case of homogeneous Neumann boundary conditions. The space X is reflexive as it is an intersection of reflexive spaces. We define where n implicitly encodes that we use scalar valued neural networks with input dimension d and arbitrary fixed depth larger or equal to two. The width of all layers (except the input and output) is set to n. With γ ≡ 0 and this definition of (A n ) n∈N the requirements of Assumption (A1) are satisfied, as can be seen by well known universal approximation results, we refer to Pinkus [27].
To proceed, the continuity of E with respect to · X is clear, hence we turn to the weak lower semicontinuity. To this end, we write E in the following form: and treat E 1 and E 2 separately. The term E 1 is continuous with respect to · X and convex, hence weakly lower semicontinuous. To treat E 2 , note that we have the compact embedding This implies that a sequence that converges weakly in H 1 ( ) ∩ L 4 ( ) converges strongly in L 3 ( ) and consequently shows that the term E 2 is continuous with respect to weak convergence in X. Finally, for fixed f ∈ X * , we need to show that the sequence (F f n ) n∈N defined in (6) is equicoercive with respect to · X . To this end, it suffices to show that the functional is coercive as it holds F f n ≥ G f . Let r ∈ R be fixed and consider all u ∈ X with G f (u) ≥ r. Then we estimate where we used the estimate due to Hölder's inequality. This clearly implies a bound on the set Description of the experiment Figure 1 shows two exemplary numerical realizations of the Deep Ritz Method with the unit disk = B 1 (0) as a domain and with right-hand sides for r 1 = 0.1 and r 2 = 0.4 corresponding to the left and right picture, respectively. Further, we considered ε = 0.01 and used a fully connected network with tanh activation and three hidden layers of width 16. Note that by Theorem 2 ReLU networks of depth log 2 (2 + 1) + 1 = 3 satisfy the universal approximation property with exact zero boundary values. Hence, the number of trainable parameters is 609 in this case. As we were solving a homogeneous Neumann boundary value problem, no penalization was needed. For the discretization of the integral over the unit disk B 1 (0), we used an evenly spaced grid and gave equal weights in the numerical approximation of the integrals to the function values at every grid point. For the optimization of the networks parameters, we used Adam with full batch size until the optimization slowed down and then used L-BFGS in order to exploit the fast local convergence properties of quasi-Newton methods. Note that in the case of f 1 , a phase transition around the ball B r 1 (0, 1/2) is energetically more favorable than the configuration in the right figure, where the radius r 2 is much larger. Remark 14 (Stability under compact perturbations) With a similar -even simpler -approach, we may also show that energies of the form fall in the Setting 5 provided E does and F is bounded from below and continuous with respect to weak convergence in X. Note also that in the space dimension d = 2 this includes the above example, however, the slightly more involved proof presented here works independently of the space dimension d.

The p-Laplacian
As an example for the uniform convergence of the Deep Ritz method, we discuss the p-Laplacian. To this end, consider the p-Dirichlet energy for p ∈ (1, ∞) given by Note that for p = 2 the associated Euler-Lagrange equation -the p-Laplace equation -is nonlinear. In strong formulation it is given by see, for example, Struwe [31] or Růžička [29]. Choosing the ReLU activation function, the abstract setting is applicable as we will describe now. For the Banach spaces, we choose where the norms · X and · B are chosen to be the natural ones. Clearly, W 1,p ( ) endowed with the norm · W 1,p ( ) is reflexive by our assumption p ∈ (1, ∞). Note that which is also reflexive. We set γ = tr, i.e., We use the same ansatz sets (A n ) n∈N as in the previous example, hence Assumption (A1) holds. Rellich's theorem provides the complete continuity of the embedding W 1,p ( ), · W 1,p ( ) → W 1,p ( ), · L p ( ) which shows Assumption (A4). As for Assumption (A3), Friedrich's inequality provides the assumptions of Lemma 8. Furthermore, E is continuous with respect to · W 1,p ( ) and convex, hence also weakly lower semicontinuous. By Poincaré's and Young's inequalities, we find for all u ∈ W 1,p 0 ( ) that Hence, a minimizing sequence in W 1,p 0 ( ) for F f is bounded and as F f is strictly convex on W 1,p 0 ( ) it possesses a unique minimizer. Finally, to provide the demicontinuity, we must consider the operator S : W Hence, the solution map S is precisely the inverse of the mapping and this map is demicontinuous, see, for example, Růžička [29].
Description of the experiment Figure 2 shows two numerical realizations of the Deep Ritz Method for the p-Laplacian with right-hand side f ≡ 1 and p 1 = 3/2 in the left picture and p 2 = 10 in the right picture. The penalization value is set to λ = 250 in both simulations to approximately enforce zero boundary values. We used fully connected feed-forward networks with three hidden layers of width 16 and GELU activation (Hendrycks and Gimpel [16]) for the left plot and ReLU activation for the right plot. The quadrature follows the same strategy as in the previous example. Note that the exact solution to the homogeneous p-Laplace problem on the disk with f ≡ 1 is given by for a suitable constant C that depends on the spatial dimension and the value of p. We see that the solution u p converges pointwise to zero for p 0 and for p ∞ the function u p tends to x → C(1 -|x|). This asymptotic behavior is clearly visible in our simulations. In case of the ReLU ansatz function, the architecture considered agrees with the analysis presented in the previous paragraph. For p 1 = 3/2, we found the GELU activation function to provide good performance. However, for the GELU activation function, establishing condition (A1) is not entirely obvious. The GELU activation is defined as GELU(x) := x (x), where is the cumulative distribution function of the Gaussian normal. It is often interpreted as a smoothed version of the ReLU since t -1 GELU(tx) → ReLU(x) and ∂ x (t -1 GELU(tx)) → ∂ x ReLU(x) pointwise for t → ∞. This gives some intuition why the GELU activation function could admit a universal approximation result with almost zero boundary values and hence satisfy (A1), however, we leave a rigorous statement for future research. This aligns very well with our numerical experiments, which do not indicate problems in resolving the zero boundary values in this practical example.

Appendix: Universal approximation with zero boundary values
Here we prove the universal approximation result which we stated as Theorem 2 in the main text. Our proof uses that every continuous, piecewise-linear function can be represented by a neural network with ReLU activation function and then shows how to approximate Sobolev functions with zero boundary conditions by such functions. The precise definition of a piecewise linear function is the following.
Definition 15 (Continuous piecewise-linear function) We say that a function f : R d → R is continuous piecewise linear, or in short piecewise linear, if there exists a finite set of closed polyhedra whose union is R d , and f is affine linear over each polyhedron. Note every piecewise linear functions is continuous by definition since the polyhedra are closed and cover the whole space R d , and affine functions are continuous.
Theorem 16 (Universal expression) Every ReLU neural network function u θ : R d → R is a piecewise-linear function. Conversely, every piecewise-linear function f : R d → R can be expressed by an ReLU network of depth at most log 2 (d + 1) + 1.
For the proof of this statement, we refer to Arora et al. [2]. We turn now to the approximation capabilities of piecewise linear functions.
Proof In the following we will denote by · ∞ the uniform norm on R d . To show the assertion, choose a triangulation T of R d of width δ = δ(ε) > 0, consisting of rotations and translations of one nondegenerate simplex K . We choose s ε to agree with ϕ on all vertices of elements in T . Since ϕ is compactly supported, it is uniformly continuous, and hence it is clear that ϕs ε ∞ < ε if δ is chosen small enough.
To show convergence of the gradients, we show that also ∇ϕ -∇s ε ∞ < ε which will be shown on one element K ∈ T and as the estimate is independent of K is understood to hold on all of R d . So let K ∈ T be given and denote its vertices by x 1 , . . . , x d+1 . We set v i = x i+1x 1 , i = 1, . . . , d to be the vectors spanning K . By the one-dimensional mean value theorem, we find ξ i on the line segment joining x 1 and x i such that Note that ∂ v i s ε is constant on all of K where it is defined. Now for arbitrary x ∈ K , we compute with setting w = d i=1 α i v i for w ∈ R d with |w| ≤ 1. Note that the α i are bounded uniformly in w, where we use that all elements are the same up to rotations and translations, where again ( * ) is uniformly small due to the uniform continuity of ∇ϕ. Noting that the W 1,∞ -case implies the claim for all p ∈ [1, ∞) finishes the proof.
We turn to the proof of Theorem 2 which we state again for the convenience of the reader.